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ABSTRACT 

Kalman  filtering  techniques  are  applied  to  a  two  sensor  bearings  only  passive 
target  motion  analysis  problem.  An  algorithm  is  developed  to  simulate  tracking  long 
range  maneuvering  airborne  targets.  The  target  tracking  performance  of  the  filter  is 
evaluated  using  computer  generated  noisy  bearing  measurements.  The  performance  of 
the  filter  is  satisfactory  given  reasonable  initial  conditions  and  measurement  noise. 


THESIS  DISCLAIMER 

The  reader  is  cautioned  that  computer  programs  developed  in  this  research  may 
not  have  been  exercised  for  all  cases  of  interest.  While  ever>'  effort  has  been  made, 
within  the  time  available,  to  ensure  that  the  programs  are  free  of  computational  and 
logic  errors,  they  cannot  be  considered  validated.  Any  application  of  these  programs 
without  additional  verification  is  at  the  risk  of  the  user. 
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I.  INTRODUCTION 

Air  defense  of  a  carrier  battle  group  is  becoming  significantly  more  complex  due 
not  only  to  the  increased  speed  and  range  of  potentially  hostile  aircraft  but  also  to 
more  capable  enemy  targeting  systems  and  greater  cruise  missile  ranges.  To  reduce  the 
probability  of  an  aircraft  carrier  being  successfully  targeted  by  an  enemy  cruise  missile 
carr\'ing  aircraft,  it  is  imperative  that  fighter  intercept  be  accomplished  beyond  the 
maximum  range  of  the  cruise  missile.  Long  range  over-the-horizon  (OTH)  target 
detection  and  tracking  are  necessary  to  achieve  this  goal. 

A  major  obstacle  common  to  all  air  defense  scenarios  is  the  enemy's  use  of 
electronic  countcrmeasures  (ECM).  Attacking  enemy  aircraft  will  undoubtedly  employ 
jarrjning  as  well  as  other  forms  of  ECM  to  degrade  or  deny  effective  tracking  by  active 
systems.  Therefore,  the  ability  to  passively  track  is  required  in  order  to  successfully 
engage  attacking  aircraft  in  a  dense  ECM  environment. 

One  viable  approach  to  this  problem  is  passive  Target  Motion  Analysis  (TMA). 
The  purpose  of  TMA  is  to  determine  the  target's  position,  course  and  speed  through  a 
series  of  passive  noisy  measurements.  For  the  air  defense  scenario,  these  passive 
measurements  may  be  lines  of  bearing  (LOB)  obtained  from  the  enemy  aircraft's 
jamming  strobes  or  from  the  electromagnetic  radiation  of  the  aircraft's  long  range 
targeting  radar. 

Passive  bearings  only  TMA  may  be  performed  by  one  or  more  sensors.  The  two 
primarv"  considerations  in  evaluating  TMA  performance  are  solution  accuracy  and 
timeliness.  Single  sensor  TMA  requires  that  the  observer  aircraft  perform  zig  zag 
maneuvers  to  establish  a  target  bearing  rate  so  that  the  range  to  the  target  may  be 
estimated.  One  drawback  to  single  sensor  T.MA  is  the  fact  that  these  maneuvers  may 
detract  from  the  observer  aircraft's  priman.'  mission.  Also,  a  reasonable  initial  estimate 
of  tlie  target's  state  (posicion,  i.uurse  aiiu  bpeeu)  is  neccsbary  lo  eiibure  liiut.  liic  irdcking 
solution  converges  in  a  timely  manner,  if  indeed  it  converges  at  all.  An  inherent 
difficulty  with  bearings  only  TMA  by  a  single  sensor  is  that  the  solution  accuracy  and 
timeliness  rely  quite  heavily  upon  a  "good"  a  priori  estimate  of  target  range. 
Consequently,  in  a  long  range  tracking  scenario  where  the  range  to  the  target  may 
exceed  several  hundred  miles,  accurate  tracking  by  a  single  sensor  using  only  bearing 
observations  is  extremely  arduous  and  rather  impractical. 


A  practical  solution  to  long  range  OTH  passive  tracking  is  multi-sensoji 
triangulation.  High  speed  air  targets  can  be  accurately  tracked  by  two  or  more  highly 
directional  sensors  that  are  spaced  sufficiently  far  apart.  The  primary  reason  that  multi- 
sensor  tracking  is  superior  to  single  sensor  TMA  is  that  estimates  of  target  range  are 
continually  being  generated  through  triangulation  of  sensor  bearing  lines.  Multi-sensor 
tracking  is  far  less  dependent  on  accurate  a  priori  state  estimates  than  is  single  sensor 
tracking  for  timely  convergence.  The  major  obstacle,  however,  in  using  the  multi-sensor 
triangulation  method  is  a  practical  one:  very  close  cooperation  is  required  between  the 
observers  in  order  to  achieve  an  accurate  tracking  solution.  Three  ingredients  are 
required  to  localize  a  target:  the  position  of  each  sensor,  the  time  of  the  observation, 
and  the  bearing  measurement  from  each  sensor.  Ideally,  all  observations  would  be 
performed  synchronously.  If  asynchronous  lines  of  bearing  are  encountered,  then 
computer  processing  is  required  to  interpolate  these  LOB's  to  produce  "synchronous" 
measurements.  The  observers  are,  in  effect,  remote  sensors  that  transmit  noisy  bearing 
data  to  a  central  processing  platform  where  the  actual  target  tracking  is  performed.  For 
trackmg  a  long  range  and  rapidly  closing  air  target,  triangulation  provides  a 
significantly  more  accurate  and  timely  tracking  solution. 

Because  each  sensor  generates  its  own  sequence  of  noisy  bearing  observations, 
the  Kalman  filter  is  ideally  suited  for  determining  a  target's  position  and  motion.  This 
thesis  investigates  the  two  sensor  bearings  only  tracking  problem  in  a  computer 
simulation  that  employs  Kalman  filtering  techniques.  The  simulation  generates  the 
target  and  observer  tracks  as  well  as  noisy  bearing  measurements  from  each  sensor  to 
the  target.  The  noisy  bearmgs  are  then  processed  by  the  Kalman  filter  to  provide 
continual  estimates  of  the  target's  state.  The  tracking  algorithm  is  used  in  several 
scenarios  to  determine  the  effect  various  sensor  bearing  accuracies  and  initial  estimates 
have  on  the  filter's  performance. 

The  aim  of  this  research  is  to  examine  how  well  the  filter  performs  in  tracking  a 
non-maneuvermg  target  before  investigating  the  more  difficult  case  of  a  maneuvering 
Largei.  The  problem  geomeiry  wiii  be  prcsenieu  firsi,  foiioweJ  by  ihe  deveiopmcni  of 
the  system  and  measurement  models.  Relevant  equations  from  Kalman  filtering  theory' 
will  then  be  briefly  reviewed  before  the  actual  tracking  algorithm  is  analyzed  in  detail. 
The  results  of  several  simulation  runs  using  various  parameters  will  be  examined.  Also, 
the  effect  of  target  maneuvers  on  filter  stability  will  be  assessed.  The  fmal  chapter  will 
summarize  the  results  of  this  research  and  will  present  conclusions  and 
recommendations  for  further  study. 
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II.  PROBLEM  DESCRIPTION 

A.       INTRODUCTION 

As  the  name  implies,  over-the-horizon  detection  and  tracking  means  positioning 
sensors  out  near  the  radar  horizon  to  look  over  the  edge  and  pass  their  observations 
back  to  a  central  data  fusion  point  for  analysis.  This  data  fusion  point  need  not  be  a 
surface  combatant;  it  could  be  a  command  and  control  aircraft.  The  use  of  an  airborne 
command  and  control  platform  extends  the  range  to  which  an  air  target  can  be 
effectively  tracked.  The  basic  idea  behind  OTH  tracking  is  to  place  the  remote  sensors 
far  enough  apart  so  that  effective  triangulation  fixes  may  be  taken  but  not  so  far  apart 
that  they  are  beyond  the  range  at  which  they  can  communicate  with  the  central 
processing  platform.  The  command  and  control  aircraft  should  ideally  be  positioned  on 
the  threat  axis  between  the  incoming  air  attack  and  the  high  value  unit  (HVU)  that  is 
being  protected.  Figure  2.1  depicts  the  general  geometry'  of  a  basic  OTH  detection  and 
trackins  scenario. 


o    \ 

\ 
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\ 

O O    _^T>.jO£ ^ 

HVU  C*C  Ale  \  IN80  £AJE*^y  A/c 

/  OF    HVU 

O     / 


Figure  2.1     Basic  Over-the-Horizon  Detection  and  Tracking  Scenario. 

In  this  chapter  the  geometn.'  of  the  two  sensor  TMA  problem  will  be  presented 
along  with  a  development  of  the  target  motion  and  noise-free  measurement  equations. 
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B.       PROBLEM  GEOMETRY 

Consider  the  target-observer  geometry  in  the  two  dimensional  plane  as  shown  in 


y^ 

i 

yr 

___jt:^ 

^' 

v^^ 

:^  T-AWET    A/C 

^i-^"""""^^'^        X^^     / 

X^  / 

/^        V/ 

1 

^"''"   ^tef«ewc*      *•                   ^*t%c9.\ 

"t 

X 

rosiTton 

Figure  2.2    Target-Observer  Geometr\'. 

Figure  2.2.  The  target  is  located  at  {xj  ,  Vj)  from  a  defined  reference  position.  The 
origin  may  be  defined  as  either  a  fi.xed  latitude/longitude  coordinate  or  the  position  of 
a  high  value  unit  such  as  an  aircraft  carrier  (whose  position  is  relatively  stationar>')- 
The  X  and  y  components  of  target  velocity  are  denoted  as  Xj  and  y-j-  and  are  the 
Cartesian  equivalents  of  the  target's  course  and  speed,  C  j  and  Vj.  Sensor  I  is  located 
at  (Xj  .  0)  and  is  only  able  to  move  along  the  x-axis  with  velocity  Xj.  Likewise,  sensor  2 
is  located  at  (0  ,  y,)  and  is  only  able  to  move  along  the  y-axis  with  velocity  y,.  As 
'>ho'''''r'  i^.  Fi^iif?  2  2  the  bep.rinc!  ^Voti  <;ijp<;of  ]  tr>  x\\p  t?!r<2c^  i^  ft  aid  th?  bearing'  frorn 
sensor  2  to  the  target  is  G^-  The  ranges  to  the  target  from  sensors  1  and  2  are  denoted 
as  r^  and  r,  respectively,  and  the  range  of  the  target  from  the  origin  is  denoted  Rj. 
I.  Problem  Assumptions 

The  following  assumptions  are  made  concerning  the  problem: 

1.  The  target  is  initially  inbound  and  remains  within  the  first  quadrant. 

2.  The  target  maintains  a  constant  speed  but  is  free  to  change  course. 
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3.  Both  sensor  positions  are  known  precisely. 

4.  Bearing  data  from  each  sensor  are  continuously  observed  and  are  received 
synchronously. 

5.  Bearing  noise  is  zero  mean  and  Gaussian  with  variances  <S-^  and  c^   for  sensors 
1  and  2,  respectively. 

6.  Target  turns  are  modeled  as  instantaneous  (i.e.,  no  turn  radius). 

7.  Target  and  sensor  altitudes  have  negligible  effect  at  long  ranges. 

The  last  assumption  is  entirely  reasonable  since  the  difference  between  the  target's 
slant  range  and  two  dimensional  range  is  less  than  a  fraction  of  1%  for  ranges 
exceeding  300  nautical  miles. 

2.  Practical  Geometrical  Considerations 

Placing  the  airborne  sensors  on  orthogonal  axes  is  chosen  not  only  because  it 
simplifies  the  geometn-'  but  also  because  it  provides  adequate  sensor  separation  with 
which  to  perform  accurate  triangulation.  Ideally,  the  most  accurate  triangulation  fix  is 
formed  from  the  intersection  of  two  perpendicular  lines  of  bearing.  Perpendicular 
LOBs  in  real  world  scenarios,  however,  are  extremely  difficult  to  obtain  for  a  number 
of  reasons.  One  reason  is  that  the  maximum  range  and  on  station  time  for  airborne 
sensors  are  limited.  .-Xlso,  if  electromagnetic  energy  from  the  target  is  being  used  to 
obtain  LOBs,  it  is  important  that  both  sensors  be  positioned  within  the  main  sector 
beam  pattern.  Figure  2.3  illustrates  some  of  these  considerations.  In  Figure  2.3  (a),  it 
can  be  seen  that  in  the  attempt  to  obtain  perpendicular  LOB's  to  the  target,  sensor  1  is 
beyond  the  radar  horizon  of  the  command  and  control  aircraft  and  is  thus  unable  to 
pass  any  bearing  observations.  Figure  2.3  (b)  shows  both  sensors  lying  within  the 
sector  scan  limits  of  the  target's  surveillance  radar.  While  it  is  not  necessar\'  for  both 
sensors  to  be  within  the  main  beam  simultaneously,  both  must  be  able  to  detect  the 
beam's  presence  vathin  a  reasonable  time  period,  a  factor  which  depends  on  the  radar's 
scan  rate. 

The  scenano  where  the  sensors  are  positioned  on  orthogonal  axes  could  easily 
be  m.ocified  to  the  more  general  case  where  the  sensors  are  located  on  radials  that  are 
separated,  say,  by  60  degrees.  Since  the  sensors  are  now  closer  together,  range 
estimates  to  the  target  would  be  somewhat  degraded.  Also,  by  adding  a  third  sensor  on 
a  radial  120  degrees  from  the  first  sensor  and  60  degrees  from  the  second  would 
provide  a  wider  sector  coverage  as  well  as  improved  tracking  accuracy. 

The  simulations  that  have  been  run  in  this  thesis  involve  extreme  ranges  from 
each  sensor  to  the  target.  It  has  been  assumed  throughout  that  the  target  and  sensor 
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Figure  2.3     Pracrical  Geometric  Consideraiions. 

aircraft  are  flying  at  medium  to  high  altitudes  so  that  all  obsen'ations  will  meet  radar 
horizon  range  constraints.  It  should  be  noted,  however,  that  the  practical  limiting 
factors  for  maximum  detection  range  are  the  strength  and  radio  frequency  (RF)  of  the 
intercepted  source  signal.  Also,  bearing  accuracies  depend  on  the  RF  of  the  signal. 
Extreme  detection  ranges,  sometimes  between  400  and  500  nautical  miles,  have  been 
used  to  represent  a  worst  case  tracking  problem;  shorter  ranges  would  yield  a  more 
accurate  tracking  solution. 

C.       SYSTEM  MODEL 

As  shown  in  Figure  2.2,  lines  of  bearing  from  two  airborne  sensors  are  used  to 
determine  the  target's  state  {position,  course  and  speed).  Using  a  Cartesian  coordinate 
system,  a  four  dimensional  state  vector,  Xj,  is  chosen. 


Xx   = 


if  J 


(eqn2.1) 


It  should  be  noted  that  the  system  model  is  in  no  way  limited  to  a  Canesian  reference 
frame  or  state  vector;  the  Cartesian  coordinate  system  was  chosen  merely  for  its 
mathematical  simplicity. 
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1.  Target  Maneuvers 

Two  basic  scenarios  are  addressed  in  this  thesis.  The  first  one  involves 
tracking  a  non-maneuvering  target  and  the  second  involves  a  maneuvering  target.  In 
both  cases  the  target  is  assumed  to  be  initially  inbound  and  any  target  maneuver  will 
consist  of  the  target  changing  only  its  course  and  maintaining  its  speed.  Figure  2.4 
shows  the  target  tracks  that  uill  be  examined  in  subsequent  chapters. 
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Figure  2.4    Representative  Target  Tracks. 

It  IS  assumed  that  target  maneuvers  can  be  modeled  by  using  white  random 
forcing  functions.  As  shown  in  Figure  2.5,  target  maneuvers  may  be  thought  of  as 
acceleration  along  its  course  (radial  acceleration)  and  acceleration  perpendicular  to  its 
course  (turn  rate).  Let  the  random  variables  6.  and  6^  denote  the  target's  acceleration 
along  its  course  and  acceleration  perpendicular  to  its  course,  respectively.  Both  6.  and 
6^  denote  random  changes  of  the  target  and  are  assumed  to  be  independent  and  zero 
mean  with  variances  c.^  and  <Ta".  Because  of  the  extremelv  lone  ranges  involved  in  the 
simulations,  target  maneuvers  have  been  modeled  as  instantaneous  changes  of  position 
according  to  the  time  interval  used.  That  is,  the  simulation  enables  the  target  to  turn 
90*^  in  two  seconds.  While  it  is  acknowledged  that  this  kind  of  turn  rate  is  quite 
artificial,  it  is  informative  to  see  what  effect  such  a  drastic  turn  rate  has  on  tracking 
performance  and  stability.  The  variances  used  in  subsequent  scenarios  are: 
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(S?  =  (300knots/sec)^ 


(eqn  2.2) 


,2  _ 


*^8    =  (45  deg/sec)' 


(eqn  2.3) 
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Figure  2.5    Geometry  of  Target  Maneuvers. 

a.    Equations  of  Motion 

Let  T  represent  the  time  interval  between  observations.  If  k  represents  the 
k^  observation  and  tj^  the  discrete  time  of  the  k^  observation,  then  T  may  be 
expressed  as 


^  ='k-  ^k-l 


(eqn  2.4) 


Referring  to   Figure  2.2   and   Figure   2.5,  target  motion  may   be  described  by  the 
difference  equations: 


i<>0  = 


'xj(k+l)" 

xj(k+l) 
yj{k+l) 

= 

yT(i^-^n 

"Xy(k)  +  T  Xj(k)  +   fi(6^,  69.  k) 

xj(k)  +  f2(6^.,  69,  k) 

yj(k)  +  T  y-j-(k)  +  f3(6^.,  6q.  k) 

yji^)  +  r4(5,.  69,  k) 


(eqn  2.5) 
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The  random  forcing  functions  fj  through  f^  are  included  to  account  for  random 
changes  in  speed  and  heading  which  occur  for  a  moving  target.  Equation  2.5  may  be 
uTitten  in  matrix  form  as 


x(k+l)  = 


or  more  concisely  as 


1    T  0    0 

~xja) 

0    10    0 

Xjik) 

0    0    1    T 

yjii^) 

0   0   0    1 

yji^) 

1^/2  0 

T  0 

0  T^I2 

0  T 


Xjik) 


(eqn  2.6) 


1  =  *^k  ^k    +    ^k^k 


(eqn  2.7) 


where 


\^  is  the  4  X  1  state  vector 


<t>^  is  the  4  X  4state  transition  matrix 

>v^  is  the  2  x  1  vector  of  random  forcing  functions 

Fj^  is  the  4x2  state  forcing  matrix. 

The  terms  of  the  random  forcing  function  vector  w^  represent  the 
accelerations  m  the  x  and  y  directions  caused  by  target  maneuvers.  The  state  forcing 
matrix  Fj.  represents  a  uniform  constant  acceleration  model  of  target  motion.  If  the 
time  interval  T  between  measurements  is  assumed  constant,  then  <J>^  may  be  replaced 
by  a  constant  state  transition  matrix  O  and  Fj^  may  be  replaced  by  a  constant  state 
forcing  matrix  F.    Revising  Equation  2.7,  the  linear  system  model  can  be  expressed  as 


^k+1  =  ^^k 


Fw, 


(eqn  2.8) 


D.      NOISE-FREE  .MEASUREMENT  EQUATION 

As  illustrated  in  Figure  2.2,  the  positions  of  sensors  1  and  2  along  with  their 
respective  bearings  to  the  target.  6^  and  6,,  uniquely  define  the  target's  position  (Xj, 
"^'j).   The  target's  position  from  noise-free  bearing  observations  may  be  expressed  as 


^    = 


(y^  cos  6j  -  Xj  sin  6j)  sin  Gj 


cos(e, 


e,) 


(eqn  2.9) 
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(x,  COS  Gj  -  y,  sin  0,)  sin  9, 
7r  = —77^ -V: ■ 


cosiGj  +  62)  (eqn2.10) 


The  positions  and  speeds  for  the  airborne  sensors  may  be  chosen  arbitrarily  for 
input  into  the  tracking  algorithm.  Each  sensor's  position  is  assumed  to  be  known 
precisely  for  each  time  interval.  The  sensors  may  both  head  inbound  or  outbound  or 
they  may  go  in  alternate  directions.  Care  must  be  exercised  in  choosing  sensor 
positions  and  speeds  so  as  to  avoid  having  lines  of  bearing  that  are  collinear  (each 
sensor  is  pointing  at  the  other).  What  results  in  this  case  is  an  extremely  thin  and 
elongated  error  ellipse  which  momentarily  degrades  tracking  accuracy  at  the  moment 
that  the  lines  of  bearing  are  coincident. 

It  should  be  noted  that  using  two  sensors  eliminates  the  need  for  any  extraneous 
observer  maneuvering  as  is  the  case  for  a  single  sensor.  The  observer  aircraft  can 
basically  fiy  straight  and  level  and  collect  more  reliable  bearings  to  the  target.  Also,  the 
sensor's  position  is  known  more  precisely  since  it  is  not  decelerating  and  accelerating 
mto  and  out  of  turns. 
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III.  KALMAN  FILTERING 

A.  INTRODUCTION 

The  technique  of  Kalman  filtering  is  ideally  suited  to  the  problem  of  passive 
tracking.  The  following  sections  briefly  describe  the  theory  and  results  of  Kalman 
filtering  and  how  it  is  applied  to  the  long  range  airborne  TMA  problem.  For  a  more  in- 
depth  development  of  the  Kalman  filter,  the  reader  is  referred  to  [Refs.  1,2]. 

B.  THE  KALMAN  FILTER 

The  purpose  of  the  Kalman  filter  is  to  keep  track  of  the  state  of  a  system 
through  a  sequence  of  noisy  measurements.  This  is  accomplished  by  recursively 
updating  an  estimate  of  the  state  by  processing  a  sequence  of  noisy  observations  in 
such  a  manner  as  to  reduce  as  much  as  possible  the  effect  of  measurement  errors. 

The  Kalman  filter  is  a  predictor-corrector  type  estimator  that  propagates  an 
estimate,  x,  of  the  target  state  along  with  an  associated  covariance  matrix,  P.  which 
reflects  the  degree  of  confidence  placed  in  the  accuracy  of  the  state  estimate.  The 
Kalman  filter  is  carried  out  in  two  alternating  stages.  First,  previous  estimates  of  x 
and  P  are  extrapolated  one  time  step  ahead  based  on  the  assumed  system  dynamics; 
this  is  referred  to  as  the  Movement  Step.  These  extrapolated  values  are  then  used  to 
compute  a  set  of  optimum  weights  called  Kalman  gains.  The  gains  are  applied  to  the 
prediction  and  to  a  new  observation  in  a  Measurement  Step,  which  provides  an 
updated  estimate  of  the  state  and  its  covariance.  This  process  is  then  repeated.   [Ref  3] 

I.  Assumptions 

The  following  assumptions  are  made: 

1.  The  random  forcing  function  w^  is  zero  mean  and  uncorrelated  with  covariance 

Qk- 

2.  The  measurement  noise  v^  is  zero  mean  and  is  correlated  with  covariance  R^. 

3.  The  random  forcing  function  >vu  and  measurement  noise  Vj^  are  uncorrelated. 

4.  The  initial  state  Xg  is  a  random  variable  with  known  mean  Xq,  j  and  covariance 
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2.  Definitions 

1.  The  estimated  state  vector  after  k  observations  is  denoted  by  x^..    and  the 
predicted  state  vector  before  the  k^  observation  is  represented  by  \|i^.i- 

2.  The  state  estimation  error  vector  £j^  is  defined  as  the  difference  between  the 
estimated  state  and  the  true  state 


^k  =  \|k  -  \  (eqn3.1) 

and  the  predicted  state  estimation  vector  z^  ^^^  is  defmed  as  the  difference  between 
the  predicted  state  and  the  true  state 

^klk-i  -  \lk-l  "  \  (^^^  ^-2) 

3.  The  covariance  of  estimation  error  matrix  Pj^jj^  is  defmed  as 

Pklk=  E(c,£j}  (eqn3.3) 

and  the  predicted  covariance  of  state  error  matrix  Pi^j^.^  is  defmed  as 

Pklk-l  =  E(£y;^.iC^,\.j]  (eqn3.4) 

4.  The  state  excitation  covariance  matrix  is  given  by 

Qk  =  ^(rni^njr'^}  (eqn3.5) 

5.  The   Kalman  filter  is  an  optimal  estimator  that  minimizes  the  sum  of  the 
variances  of  the  estimation  error,  i.e., 

E{Ci(k)-}  +  E(e2(k)2}  +    .  .  .    +  E{e„(k)2}  (eqn  3.6) 
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Enter  known  matrices  and  a  priori  estimates: 
%|.l)  •  P(0|-1)  •   R(O)  .   H  ,  (P 

Compute  the  Kalman  gain: 

MEASUREMENT  STEP 

^(k|k)  =  ^(k|k-l)  -^  ^(k)'  ^(k)  -  "%|k-l)  ) 
P(k|k)  =  t  I  -  %)H}P(t^|t.i) 

MOVEMENT  STEP 


^(k+Ilk)  -  ^^(k!k) 
P(k+l|k)  =  ^P(k|k)^^   ^    Q 


(k) 


Compute  R/j^x  and  Q^-j^n 
Increment  k  bv  l 


Figure  3.1     Kalman  Filter  Algorithm. 
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3.  Kalman  Filter  Algorithm 

Figure  3.1  summarizes  the  discrete  Kalman  filter  algorithm.  For  the  particular 
T.MA  problem  presemed  in  this  thesis,  the  2  x  4  measurement  matrix  H,^  and  the  4  x  4 
state  transition  matrix  Oj^  are  both  known,  constant  matrices  and  may  be  represented 
by  H  and  O.  An  a  priori  estimate  Xq,  ^  of  the  target's  state  with  an  associated  initial 
error  covariance  matrix  Pgi.!.  as  well  as  an  initial  estimate  of  the  measurement  noise 
covariance  matrLX  Rq  must  be  input  into  the  filter  algorithm.  The  algorithm  computes 
the  Kalman  gain  Gj^  based  on  these  a  priori  values  and  then  updates  the  estimate  of 
the  target's  state  when  it  receives  a  measurement.  The  error  covariance  matrix  is  also 
updated.  Next,  the  state  estimate  and  its  error  covariance  matrix  are  projected  one  time 
step  ahead  based  on  the  assumed  system  dynamics.  The  measurement  noise  covariance 
Rj^  and  the  state  excitation  covariance  matrix  Q,^  are  then  computed  before  k  is 
incremented  by  one  and  the  whole  process  is  repeated. 

C.       FUNCTIONS,  MATRICES,  AND  EQUATIONS 

In  this  section,  the  Kalman  filter  algorithm  will  be  applied  to  the  long  range 
passive  airborne  tracking  problem.  A  brief  derivation  of  the  random  forcing  function 
>V|^,  the  state  excitation  covariance  matrix  Q^,  the  measurement  equation  z^.  and  the 
measurement  noise  covariance  matrix  R^  is  given  next. 

1.  Random  Forcing  Function 

Recalling  equation  (2-6),  the  two  dimensional  random  forcing  function  w^ 
represents  the  acceleration  in  the  x  and  y  directions  caused  by  target  maneuvers. 


)^k      = 


(eqn  3.7) 


where  v,    =  (  x.  ^  4-  y  ^  )l/2_ 


Since  the  random  variables  6^.  and  6^  were  assumed  to  be  zero  mean,  it 
follows  that  the  random  forcing  function  w^^  is  also  zero  mean.  The  variances  of  the  x 
and  y  accelerations,  denoted  by  G-^  and  (T---  respectively,  are 

X  y 


(eqn  3.S) 
(eqn  3.9) 
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The  covanance  of  the  x  and  y  acceleration  c--^  is 

xy 


V  =   ^[^K  yj   =    -x,y,^'  *  iLk4  (eqn  3.10) 


Therefore,  the  random  forcing  function  covariance  matrix  Q.  '  is 


Q.    =  E{»,  w/  } 


T  X  _ 


k     k 


/  ^ 


« 


^.  A 


»*» 


(eqn  3.11) 


where  (T-.'  ,  (T-^  .  and  (S-.?  are  computed  at  the  predicted  values  of  Xy  and  vt-. 

A  y  av  1  "'X 

2.  State  Excitation  Covariance  Matrix 

The  purpose  of  the  state  excitation  covariance  matrix  Qj.  is  to  account  for 
model  inaccuracies  or  for  a  target  that  has  maneuvered.  It  is  basically  a  "procedure  for 
masking  the  effects  of  modeling  errors"  [Ref.  2:  p.  163].  In  effect,  the  state  excitation 
covariance  matrix  increases  the  size  of  the  predicted  covariance  of  error  matrix  which 
in  turn  increases  the  filter  gains.  As  more  observations  are  processed,  Qj^  prevents  the 
Kalman  gains  from  approaching  zero  by  continually  injecting  uncertainty  into  the 
predicted  state  estimate  at  each  iteration.  A  nonzero  Qj^  slightly  degrades  the  filter's 
accuracy  when  the  target  is  not  maneuvering  but  it  helps  prevent  filter  divergence  when 
the  target  does  maneuver.  As  stated  in  equation  3.5,  the  state  excitation  covariance 
matrix  is 


Qk  =  '"Q,rT  = 


f-» 


t4  r.i 

2.     «  M     «Y 


SYmmETRiC 


i"-^. 
^-^ 


2.      "1 

TV, 


(eqn  3.12) 


3.  Measurement  Equation 

In  this  TM.A.  problem,  the  observations  are  noisy  (x,y)  positions.  It  is  the 
intersection  of  noisy  lines  of  bearing  that  form  the  noisy  (x,y)  position  of  the  target 
that  is  input  into  the  Kalman  filter  algorithm.  Because  the  observations  are  of  the 
same  form  as  the  state  vector,  the  measurement  equation  is  linear  and  is  expressed  as 
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\=^\^\ 


(eqn3.13) 


where 


Zj^  is  the  2  X  1  measurement  vector 
H  is  the  2  X  4  measurement  matrix 
Vj^  is  the  2  x  1  measurement  noise  vector 
Xj^  is  the  4  x  1  state  vector. 


The  equation  may  be  ^\'Titten  explicitly  as 


x 

1    0   0   ol 

\ 

0    0    10 

>'k 

>k 

'Ik 

''2k 


(eqn3.14) 


The  most  important  part  of  the  measurement  equation  is  an  accurate  description  of  the 
measurement  noise  vector  v^.  The  measurement  noise  vector  expresses  the  statistical 
nature  of  the  noisy  {x,y)  position  that  is  derived  from  the  intersection  of  two  noisy  lines 
of  bearing.  These  bearing  measurement  errors  are  assumed  to  be  independent  and  zero 
mean  with  variances  (T^^  and  (T2",  for  sensor  1  and  sensor  2  respectively. 

It  is  important  to  note  that  the  bearing  errors  between  sensors  are  statistically 
uncorrected;  one  sensor's  bearing  accuracy  has  nothing  to  do  with  any  other  sensor's 
bearing  accuracy.  However,  in  describing  the  resulting  intersection  in  Cartesian 
coordinates,  the  noisy  x  and  noisy  y  positions  are  correlated.  The  only  case  where  the 
noisy  x  and  noisy  y  positions  are  uncorrected  is  when  the  lines  of  bearing  are 
perpendicular. 

4.  Error  Ellipses 

An  intuitive  way  to  visualize  the  measurement  equation  is  through  the  concept 
of  error  ellipses.  Error  ellipses  give  a  geometric  picture  of  the  region  around  a  noisy 
position  or  estimate  where  the  true  value  is  considered  to  lie.  Figure  3.2  shows  a 
Hiv^rirtte  G3u*^sinn  rrobcibili'"^'  dcn'^if'  function  f^rm?d  b'''  th?  !ntcr'^?ct!cr.  cf  f^'o  I'.nes 
of  bearing  with  independent  Gaussian  distributions. 

As  can  be  seen  in  Figure  3.2,  the  lines  of  bearing  intersect  at  an  oblique  angle, 
forming  an  asymmetric  hump.  While  the  bivariate  Gaussian  probability  density 
function  gives  an  interesting  three  dimensional  depiction  of  two  normally  distributed 
bearing  errors,  it  does  not  provide  the  information  that  is  really  needed,  quickly.  What 
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Bivariate  Normal  Density  Surface:  Sx  =  3,  Sy=5.  p=.5 


A^0 


^^^      X  AXlS°°n  NWt) 


Figure  3.2     Bivariate  Gaussian  Probability  Density  Function. 

is  needed  is  an  accurate  picture  of  the  measurement  (or  estimation)  error.  This 
uncertainty  is  best  expressed  geometrically  by  the  error  ellipse.  The  term  "error  ellipse" 
refers  to  the  two  dimensional  surface  of  constant  probability  density.  Figure  3.3 
presents  these  error  ellipses  as  contour  hnes  of  the  bivariate  Gaussian  probability 
density  function  shown  in  Figure  3.2. 

The  various  ellipse  sizes  in  Figure  3.3  correspond  to  different  constant 
probabilities.  The  fact  that  the  ellipses  are  also  rotated  imphes  that  the  uncertainty  in 
measurement  error  is  indeed  correlated  with  respect  to  x  and  y.  The  actual  probabilities 
within  a  specified  error  ellipse  may  be  computed  through  lengthy  integration  of  the 
bivariate  Gaussian  probability  density  function  over  the  surface  of  the  ellipse.  Some 
computed  probabilities  of  the  true  value  lying  within  the  Iff.  2cr,  and  3<T  error  ellipses 
are  .394,  .865.  and  .989  respectively  [Ref.  4:  pp.  4-49]. 

Error  ellipses  are  extremely  useful  in  examining  postion  error.  Matrices 
containing  the  x  and  y  position  terms  convey  analytically  what  error  ellipses  display 
graphically.  A  2  >=  2  error  covariance  matrix  which  contains  position  components  x  and 
y  is  able  to  completely  describe  an  ellipse.  The  main  diagonal  terms  represent  the 
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Contours  of  the  Bivoriote  Normal  Density  Surface 
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Figure  3.3     Error  Ellipses  as  Contour  Lines. 

variances  of  the  x  and  y  positions  respectively.  The  ofT  diagonal  terms  represent  the 
degree  of  x-y  coupling  and  the  orientation  of  the  error  ellipse  in  the  x-y  plane. 
5.  Covariance  of  Measurement  Error  Matrix 

The  covariance  of  measurement  error  matrix  Rj^  uses  the  concept  of  error 
ellipses  to  accurately  describe  the  noisiness  and  degree  of  coupling  of  (x,y) 
measurements  obtained  from  intersections  of  noisy  LOB's.  The  terms  of  the  covanance 
of  measurement  error  matrix  Rj^  depend  on  the  standard  deviations  of  bearing  error  <ij 
and  Cj  of  sensors  1  and  2  as  well  as  the  angle  at  which  the  lines  of  bearing  intercept. 
The  covariance  of  measurement  error  may  be  expressed  as 


R  = 


g;.  sin 9^  •!■  q^-  e-oi9,  >q:^  t.^e,  coiB^  -  a-J cosO.  <.^g, 


coi-(e,-»6f) 


-ff",  Smdxt^idf^    —  OJ.   <.os&,  Stnd, 


coi 


■(^.^e^) 


cos 


■(e.'>B^) 


(;.n3-/5) 
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where  0j  and  9^  denote  noisy  bearing  observations  from  sensors  1  and  2,  respectively. 
The  subscript  k  has  been  intentionally  deleted  from  equation  (3-15)  only  for  ease  of 
notation.  At  each  discrete  time  interval  tj^,  new  values  of  0j  and  62  are  generated  "with 
which  to  compute  the  new  measurement  error  covariance  matrix,  R^  A  complete 
derivation  of  equation  (3-15)  is  given  in  Appendix  A  for  the  interested  reader. 
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IV.  THE  ALGORITHM 

A.       INTRODUCTION 

This  section  discusses  the  development  of  the  tracking  algorithm.  The  algorithm 
is  designed  to  simulate  tracking  a  long  range  inbound  enemy  air  target  by  triangulating 
noisy  bearing  observations  from  two  airborne  sensors.  We  want  to  be  able  to  track  a 
non-maneuvering  target  within  a  one  percent  range  error.  For  a  maneuvering  target, 
we  desire  a  stable  filter  response  which  quickly  converges  to  the  target's  new  state.  The 
effect  of  various  sensor  bearing  errors,  a  priori  state  estimates  and  initial  error 
covariance  matrices  on  filter  accuracy  and  convergence  time  are  investigated.  Also,  the 
effect  of  target  maneuvers  on  filter  stability  is  analyzed.  Basically,  the  algorithm 
performs  three  functions: 

1.  The  target  and  sensor  tracks  are  generated. 

2.  Noisy  bearing  observations  are  simulated  using  a  random  number  generator. 

3.  The  noisy  measurements  are   processed  by  the   Kalman  filter  algorithm  to 
generate  estimates  of  the  target's  state. 


B.  TARGET  TRACK 

As  mentioned  m  Chapter  2,  three  target  track  scenarios  are  investigated: 
nonmaneuvermg.  gentle  turn,  and  hard  turn.  In  all  three  cases  the  initial  target  position 
IS  (410  nmi,  430  nmi)  with  X  and  V  velocities  of  -400  knots  and  -380  knots 
respectively. 

C.  GENERAL  SIMULATION  SCENARIO 

The  overall  purpose  of  this  simulation  is  to  be  able  to  track  an  inbound  enemy 
aircraft  before  it  flies  within  300  nmi  of  the  high  value  unit.  By  using  two  sensors  which 
circ  ebbeniiall)  dbie  lo  "peek"  over  the  norizon,  a  long  range  OTH  ilgnter  intercept  of 
the  target  aircraft  may  be  accomplished.  For  all  simulation  runs,  an  initial  target  range 
of  600  nmi  is  chosen.  This  allows  the  sensors  to  passively  track  for  thirty  minutes,  and 
it  enables  fighter  intercept  to  occur  beyond  450  nmi  of  the  HVU,  depending  on  the 
fighter's  initial  position  and  fuel  state.  The  target  aircraft  is  also  assumed  to  be  flying 
inbound  at  approximately  600  knots. 
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1.  Algorithm  Flow 

The  algorithm  can  be  broken  down  into  the  following  steps: 

1.  Define  the  true  target  track. 

2.  Define  the  observer  tracks. 

3.  Enter  a  priori  estimates  XQip   PqI-I'   ^^^  ^^-    Enter  bearing  error  variances  <Tj^ 
and  Cj. 

4.  Compute  the  noise  free  bearings  from  each  sensor  to  the  target  for  each  time 
interval. 

5.  Compute   random   sensor   bearing   errors   using   computer   generated  normal 
distribution. 

6.  Add   the  random  bearing  errors   to   the  noise   free  bearings   to   create  noisy 
bearing  measurements. 

7.  Compute  the  noisy  (x.y)  position  that  results  from  the  intersection  of  two  noisy 
lob's. 

S.      Input  this  noisy  (x,y)  measurement  into  the  Kalman  filter  algorithm. 

D.  T.ARGET  TRACK 

As  mentioned  in  Chapter  2,  three  target  track  scenarios  are  investigated: 
nonmaneuvering,  gentle  turn,  and  hard  turn.  For  all  scenarios,  the  initial  target 
position  is  (410  nmi,  430  nmi)  with  x  and  y  velocities  of  —400  knots  and  -  380  knots 
respectively.  For  the  scenarios  where  the  target  maneuvers,  a  value  is  input  for  the  time 
chat  the  maneuver  is  to  take  place.  .Also,  values  for  the  x  and  y  velocities  are  input  for 
the  second  leg  of  the  target  track.  It  should  be  noted  that  the  turn  radius  of  the  target 
is  no:  taken  into  account  in  the  target  track,  for  at  the  extreme  distances  being 
investigated,  target  maneuvers  almost  appear  as  point  turns. 

E.  OBSERVER  TRACKS 

The  observer  tracks  each  begin  at  295  nmi  on  their  respective  axis  and  travel 
inbound  at  420  knots.  These  observer  tracks  were  chosen  to  be  inbound  to  represent  a 
more  realistic,  worst  case  type  scenario.    Ideally,  it  is  desired  to  have  both  observer 

constraints   and   maximum   on-station   time   for  the   observer  aircraft   are  important 
practical  considerations  that  must  be  taken  into  account. 

F.  INITIALIZATION 

In  the  first  senes  of  simulations,  different  combinations  of  the  initial  state 
estimate  .\qi  j  and  the  initial  error  covariance  matrix  Pqi  ^  are  tested.  For  the  first 
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simulation,  and  the  one  by  which  the  other  simulations  are  compared,  the  a  priori  state 
estimate  is 


2io|.| 


400  nmi 
420  knots 
400  nmi 
420  knots 


and  its  associated  initial  error  covariance  matrix  is 


(50nmi)2 

0 

0 

0 

0 

(50kts)- 

0 

0 

0 

0 

(50nrrd)' 

0 

0 

0 

0 

(50kts 

P.H  = 


The  initial  measurement  error  covariance  matrix  for  the  one  degree  bearing  error  case 
is 


R.  = 


(7nmir 
(5nmi'r 


(5nmi)^ 
(7nmi)' 


G.       NOISY  BEARING  GENERATION 

The  Box-Miiller  method  is  used  to  generate  normally  distributed  bearing  errors. 
Basically,  it  is  a  mapping  technique  that  uses  an  algebraic  identity  to  establish  a  one  to 
one  relationship  between  a  uniform  random  variable  and  a  normal  random  variable. 
Two  random  L'(0,1)  numbers,  Uj  and  L'^^are  transformed  into  independent  N(0,1) 
random  numbers,  Nj  and  N-,  using  the  equations 

Nj  =  (-21nL"i)'"cos27rU2 

N2  =  (-21nU^)^'-  sm27tL'2 
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Figure  4.1     Normally  Distributed  Beanng  Error. 

Figure  4.1  presents  a  histogram  shov^'ing  the  normal  distribution  of  the  bearing  error. 
These  normally  distributed  random  numbers  are  then  multiplied  by  the  standard 
deviation  of  measurement  error  for  each  sensor  to  produce  two  independent  normally 
distributed  bearing  errors  for  6.  and  9^. 
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V.  SIMULATION  RESULTS 

A.  INTRODUCTION 

The  purpose  of  this  chapter  is  to  show  through  various  scenarios  the  effect  of 
different  initial  state  estimates  and  measurement  noise  levels  on  the  stability,  accuracy, 
and  convergence  time  of  the  algorithm.  In  the  following  pages,  fourteen  simulations 
involving  three  scenarios  are  presented.  As  shown  in  Figure  2.4,  the  three  scenarios 
include  a  nonmaneuvering  target  track,  a  target  track  with"gentle  turn,  and  a  track  with 
a  hard  turn.  The  first  scenario  provides  the  reference  with  which  other  cases  may  be 
compared.  Unless  otheru^ise  noted,  all  simulations  use  a  two  second  time  interval 
between  measurements.  Also,  all  of  the  simulation  results  depict  the  cases  of  one 
degree  and  five  degree  sensor  bearing  errors.  It  should  be  pointed  out  that  in  order  to 
isolate  the  effect  of  various  parameter  changes,  a  single  random  number  seed  is  used 
throughout  to  represent  a  specific  noise  histor>'.  There  has  been  no  attempt  to  create 
statistics  based  on  ensemble  of  noise  histories  due  to  the  extreme  computational  time 
required. 

B.  TYPES  OF  GRAPHS 

Five  graphs  are  used  in  all  simulations.  These  include  the  x  and  y  positional 
errors,  the  x  and  y  velocity  errors,  and  the  percent  range  error.  For  the  case  of 
positional  errors,  the  updated  state  estimate  of  x  and  y  position  is  subtracted  from  the 
target's  true  x  and  y  position.  Likewise,  for  the  case  of  velocity  errors,  the  updated 
state  estimate  of  velocity  is  subtracted  from  the  target's  true  x  and  y  velocity.  Range 
percent  error  is  computed  as 


Ir  -r  / 

Range  Percent  Error  =      ^    ^     x   100 

where  Rj  denotes  the  target's  true  range  from  the  origin  and  Rj  is  the  updated 
estimate  of  target  range  using  the  updated  state  estimate  for  the  x  and  y  positions.  In 
some  simulations,  the  measurement  residual  error  is  plotted.  The  measurement  residual 
is  defi..  d  as  the  dilTerence  between  the  actual  noisy  measurement  and  the  predicted 
state  estimate.    It  may  be  expressed  as  the  quantity 
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C.       SCENARIO  1 

Scenario  1  consists  of  ten  simulations  that  demonstrate  the  effect  of  various  a 
priori  estmiates,  measurement  noise  levels,  and  time  intervals  between  measurements 
on  filter  performance.  In  all  ten  simulations,  the  target  is  on  constant  course  and 
speed.  Figure  5.1  illustrates  the  target's  true  track  along  with  noisy  measurements. 
Note  by  the  scale  that  only  a  portion  of  the  first  quadrant  is  depicted.  In  the  first 
simulation,  the  initial  difference  between  the  true  target  state  and  the  initial  state 
estimate  is 


2<T  ~  2Soj-i     = 


10  nmi 
20  kt 
30  nmi 
40  kt 


and  the  a  priori  error  covariance  matrix  is 


■oi-i 


(5o«-..)    O        O      O 

o      o  (so.Jo 


(Sokij 


Overall,  it  can  be  seen  that  for  the  one  degree  bearing  error  case,  the  algorithm 
tracks  quite  well.  Referring  to  Figure  5.1,  the  x  velocity  error  initially  gets  worse  before 
it  gets  better.  It  is  not  until  after  five  minutes  have  elapsed  that  the  x  velocity  error  is 
less  than  the  a  priori  estimate.  As  can  be  seen  from  Figure  5.1,  the  tracking  accuracy 
for  the  five  degree  bearing  error  case  is  fair.  The  one  degree  bearing  error  case 
definitely  demonstrates  quick  and  accurate  filter  convergence;  the  five  degree  bearing 
case  seems  to  meander  almost  randomly.  The  sudden  rise  for  the  five  degree  case  seems 
to  be  an  anomalous  disturbance.  The  x  and  y  velocity  gains  are  directly  related  to 
beanng  accuracy;  the  higher  the  accuracy,  the  greater  the  gain.  For  the  five  degree 
bearing  error  case,  the  velocity  gain  never  exceeds  0.2. 

The  following  nine  cases  are  basically  variations  of  the  first  case.  Table  1  lists 
the  parameters  that  have  been  changed  for  each  case. 
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TABLE  1 
SCENARIO  PARAMETER  CHANGES 


FIGURE 

r/iRcer 

T(UCK 

rime 

riVTEIlvAL 

/\  r 

coK^>^e*4T% 

S.l*?.Z 

2-c 

W;3<?   -V20    fw?  -Yiff 

5.3 

5TRM/<«r 

5«c 

VOO    -Y2o    *«IOO  -VZO 

/NCACASE»    T/mf     //MTKCvZ/^t 

S.H 

SnWl6MT 

25« 

VJ«   -5fO    I20    -ISO 

[■■"! 

5.5 

STRAI^KT 

a^^c 

VJO    -SVO    I20    -ISO 

<«rf«e    ««  •i.vtf  ;     ^     m<fr;i,  uncorreUi^J 

S.C 

STRAIGHT 

25*< 

tfoo    -tfo   Moo  -see 

P,or  A^Or."    Ve/oe/^y,  ^-o./  p«*.-'.»n 

57 

STKKtGva 

2s^c 

Zoo  -V20  V60  -va? 

poor  «  priori    poii-^'ow,    ^ce^  vc/*c<°'/y 

5.8 

STRAIGHT 

a«. 

ioo   -««o  VW  -foo 

Poor   A    prior,    po«i'T.o»>   aj«^  Vc/oeify 

S.I 

STIWitfMT 

2s*c 

2o<3    -<W    VgO  -500 

5.^,  *,*-Ue.    P„.,=  2?ool 

SAO 

SnwitfHT 

25«C 

*r(?0     -V<0    V2«  -370 

Verw   «•««/  «  fr»'»rt'  e»f.  J          ^m'^^^ 

s.n 

snwi^UT 

2s*< 

«A)0    -VW  Y«o  -Vto 

)«..«  «f   RtS   e«»«   KHctfi   ^,-^,  *;*' 

5.12 

ficisrru  TUKH 

25«r 

Hoo   -Vto  i^  -Y2r> 

•         -/O0»tT 
V^*^    'S-vn;    new    y-.,^^ 

5./3 

HARD  TuKhJ 

Zuc 

V#o   -HZO  Ytw  -YW 

5./V 

fi€AJT|.E  TlHtiJ 

2fcc 

V«o    -ilo  Hoo  'iZo 

S<u.«  At  ;*  f,'^  S.lt,  '*i  •"•*-f»    fr*ti.'nj 
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TRUE  X  -  FILTERED  X  POSmON  VS.  TIME 
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Figure  5.1     Straight  Track  Reference  Case. 
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X  AND  Y  VELOCITY  GAINS  VS.  TIME 


X  AND  Y  POSITION  GAINS  VS.  TIME 
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Figure  5.2     Straight  Track  Reference  Case  (Cont'd). 
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TRUE  X  -  FILTERED  X  POSITION  VS.  TIME 


TRUE  Y  -  FILTERED  Y  POSmON  VS.  TIME 
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TRUE  X   -   FILTERED  X  VELOCITY  VS.   TIME 
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Fi2ure  5.3     Time  Inten'al  =  5  sec. 
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TRUE  X  -  FILTERED  X  POSITION  VS.  TIME 
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Figure  5.4    APXHAT  =  430  -540  120  -50,  with  R  Matrix. 
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TRUE  X  -  HLTERED  X  POSITION  VS.  TIME 


TRUE  Y  -  FILTERED  Y  POSITION  VS.  TIME 
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TRUE  X  -  FILTERED  X  VELOCITY  VS.  TIME 
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Figure  5.5     Same  as  Case  3,  Uncorrelated  R  Matrix. 
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TRUE  X  -  FILTERED  X  POSfTION  VS.  TIME 
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TRUE  Y  -   FILTERED  Y  VELOCITY  VS.  TIME 
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Figure  5.6    APXHAT  =  400  -150  400  -500. 
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TRUE  X  -  FILTERED  X  POSITION  VS.  TIME 
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Figure  5.7    APXHAT  =  200  -420  480  -420. 
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TRUE  X  -  FILTERED  X  POSmON  VS.  TIME 
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Figure  5.8    APXHAT  =  200 -150  480 -500    P=  10^1. 
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TRUE  X   -  FILTERED  X  POSITION  VS.  TIWE 
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Figure  5.9     APXHAT  =  200  -150  480  -500     P=  25001. 
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TRUE  X  -  FILTERED  X  POSITION  VS.  TIME 
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Fi2ure5.10    APXHAT  =  400 -410  420 -370     P=  1001. 
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Figure  5.11     (Tj  =  5  deg   (T^  =  1  deg. 
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Figure  5.12    Gentle  Turn  With  Q  Matrix. 
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TRU£  X  -   FILTERED   x   POSITION  \'S.  T;m: 
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Figure  5.13     Hard  Turn  With  Q  Matrix. 
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Figure  5.14    Case  II  For  45  min.  run. 
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VI.  CONCLUSIONS  AND  RECOMMENDATIONS 

A.       CONCLLSIONS 

The  purpose  of  this  thesis  was  to  investigate  the  two  sensor  bearings  only  passi\'e 
tracking  problem  using  Kalman  filtering  techniques.  A  computer  simulation  was 
de\eloped  to  generate  the  target  track  and  the  noisy  bearing  observations  from  each 
sensor.  Filter  performance  was  e.xceptional  for  the  nonmaneuvering  target  case.  With 
one  degree  sensor  bearing  error  and  two  second  measurement  intervals,  the  filter  was 
able  to  consistently  track  the  target  to  within  one  quarter  of  one  percent  range  error  in 
the  first  fi\e  minutes.  As  was  expected,  filter  accuracy  was  degraded  as  bearing  error 
was  increased.  The  tracker  performed  reasonably  well  for  sensor  bearing  errors  as  high 
as  eight  degrees. 

For  the  case  of  a  maneuvering  target,  filter  performance  was  marginal.  Filter 
convergence  to  an  accuracy  attained  prior  to  the  target  maneuver  did  not  occur.  The 
use  of  a  state  excitation  covariance  matrix  by  itself  was  not  sufficient  enough  to 
properly  account  for  target  maneuvers.  What  is  needed  is  a  reliable  zig  detector  that 
quickly  recognizes  target  maneuvers  so  that  the  filter  gains  may  be  reinitialized.  The 
problem  of  detecting  target  maneuvers  is  not  trivial.  At  long  ranges,  error  ellipses  may 
be  twenty  to  forty  times  larger  than  the  actual  distance  the  target  has  moved.  Sifting 
out  a  bona  fide  target  maneuver  from  these  extremely  noisy  measurements  is  quite 
difficult.  Determining  a  target  maneuver  by  attempting  to  find  a  pattern  in  the 
measurement  residuals  is  only  successful  if  the  time  between  measurements  is  greater 
than  thirty  seconds.  A  drawback  to  this  approach  is  that  filter  accuracy  is  degraded 
because  fewer  measurements  arc  being  processed. 

From  the  simulation  runs  it  was  found  that  the  most  influential  factors  m 
determining  tracking  accuracy  and  speed  of  convergence  were  the  sensor  bearing 
accuracies  and  the  time  between  measurements.  Factors  that  contributed  to  a  lesser 
extent  included  the  accuracy  of  the  initial  state  estimate  along  with  its  associated 
degree  of  confidence  and  the  positions  of  the  sensors.  Inaccurate  a  priori  information 
did  not  degrade  the  filter's  accuracy,  it  only  increased  the  time  for  convergence.  Filter 
accuracy  improved  as  the  lines  of  bearing  came  closer  to  being  perpendicular. 
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B.       RECOMMENDATIONS 

This  study  is  by  no  means  complete.  Some  areas  for  further  study  include  the 

following: 

1.  Run  an  entire  ensemble  of  simulations  (perhaps  1000  runs)  to  generate  reliable 
statistics  on  the  filter  s  performance. 

2.  Investigate  more  fully  a  method  to  detect  target  maneuvers  so  that  adaptive 
control  techniques  may  be  used  to  alleviate  the  problem  of  filter  divergence. 

3.  Examine  the  utility  of  using  more  sensors  to  cover  a  comparable  sector.  Are 
more  sensors  necessarily  beneficial.' 

4.  Peribrm  the  simulation  using  a  diilerent  coordinate  system  (such  as  polar 
coordinates)  and  compare  the  results  to  those  obtained  using  a  Cartesian 
model. 
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APPENDIX  A 
DERIVATION  OF  .MEASUREMENT  ERROR  COVARIANCE  MATRIX 

In    [his    appendix    the    measurement    error    covanancc    matrix    Rj.    is    derived. 
Recalling  the  measurement  equation  from  Chapter  3.  the  measurement  noise  vector  v 
is  assumed  to  be  Gaussian  and  zero  mean  with  variance  R.  .  That  is, 


V,  -  N(0,R,  ) 


(A-\) 


The  purpose  of  R^  is  to  statistically  describe  the  noisiness  of  the  x  and  y  measurements 
obtained  through  the  intersection  of  noisy  lines  of  bearing.  Basically,  the  2  by  2 
micasurement  error  covariance  matrix  describes  this  noisiness  by  displaying  the  variance 
and  covariance  of  the  noisy  x  and  y  measurements  in  terms  of  each  sensor's  bearing 
measurement  and  accuracy.  Referring  to  Figure  A.l,  the  position  (Xj.Yj)  represents 
a  possible  true  position  of  the  target  based  on  noisy  sensor  bearing  observations,  G^ 
and  6,.  The  position  of  the  target  (Xj,  Yj)  is  a  jointly  distributed  random  variable 
whose  expected  value  coincides  with  the  intersection  of  the  bearing  lines. 
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Figure  A.l     Target  Observer  Geometr>'  Revisited. 
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To  develop  the  relationship  between  Xj  and  Yj.  the  LOBs  from  each  sensor 
may  be  expressed  in  the  general  form  for  the  equation  of  a  line: 


y     sm  ^,     ••■    y  cos  S^      -     X,  sinP^     -     O 


X      to$l 


(A-Z) 
(A-i) 


The  distance  from  the  position  (X j,  Yj)  to  each  sensor  Ime  of  bearing  is  denoted  by 
d,  and  d^  for  each  sensor  LOB  respectively.  From  the  problem  geometn."  and  by  usmg 
the  equation  for  the  distance  between  a  pomt  and  a  line,  the  displacement  distances  d, 
and  d^  may  be  expressed  as 


(A'5) 


Since  both  sensor  bearing  errors  are  assumed  to  be  Gaussian  zero  mean  random 
variables  with  bearing  variances  (Tj^  and  (T^^,  it  follows  that  the  displacement  random 
variables  d^  and  d-,  are  also  zero  mean  Gaussian  with  displacement  variances  dj?  and 
or,/  respectively.  Figure  A. 2  illustrates  the  relationship  between  the  displacement 
variance  and  the  sensor  bearing  error  variance. 


H' 

-? 

^<»; 

^ 

X, 

S€niet  1 

X 

Figure  A. 2     Relationship  Between  Bearing  and  Displacement  Errors. 


In  this  figure,  the  bearing  error  standard  deviation  for  sensor  1,  Cy  is  expressed 
in  degrees  or  radians  whereas  the  displacement  error  standard  deviation  for  sensor  1.  a, 
,  is  expressed  in  nautical  miles.  They  are  related  by 

(Tj,  =  rj  tand,  (fy-c,) 

where  rj  is  the  approximate  distance  from  sensor  1  to  the  target  in  nautical  miles. 

Therefore,  the  displacement  distances  dj  and  d^  are  normal  random  variables  that 
may  be  v/ritten  as 

d,  ~  N(0,  <T,?)  ^^.g) 

Having  described  the  displacement  random  variables  d^  and  d^  from  sensor  1   and 
sensor  2  lines  of  bearing,  equations  A-4  and  A-5  may  be  rewritten  as 

Xr  ex,,  e,  +  Yt  s..  0,   =  yz  si.e^  ^  N^  (A-io) 

where    N^   =   N(0,  G^?)  and  N^  =   N(0.  a  J).    Solving  equations  (A-9)  and  (A- 10)  for 
X-r  and  "i'y  yields; 


or 


where 
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and 


f{    -    •"^'''^t 


^  -  COT  e, 


Note  that  since  Xy  and  Yy  are  linear  combinations  of  normally  distributed  random 
variables,  they  must  also  be  normally  distributed.   That  is. 


Xj  -^  N(  Mx^'  ^xi) 


where        ^*X  ~  "^l 


*^xi  =   ^2'  ^I'^  +   •'^3^  <^2'^ 


(A'\S) 


and  ^Y  =  B, 


ay-  =  B.-  <5^  +  83^  Q^ 


Now,  the  measurement  error  covariance  matrix  R  may  be  written  as 


(A'\C) 


R  = 


Var(Xj)  Cov(Xj,  \j) 

Cov(Xj.  Yj)     Var(Yj) 


(^-n) 


By  definition,  the  covariance  of  the  random  variables  Xj  and  \j  is  the  expected 
value  of  their  product  minus  the  product  of  their  means.   That  is, 


Cov(X-j-,  Yy)  =  e[Xj  Yj]   -  Mx^My 


(A-\8) 


Substituting  equations  (A-Il)  and  (A-12)  into  the  above  equation  yields 
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By  noting  that  Nj  and  N-,  are  zero  mean,  the  expected  value  of  their  product  is  the 
product  of  their  means.  That  is,  E[NjN\7  =  C  [Ni]  e|^->/  ^  ^-  t-smg  this  fact, 
the  covanance  between  the  random  variables  Xj  and  \j  is  simplified  to 


Cov(Xy.  Yj)  =  A,  B,  (Tj,-  +  A3  B3  (7 J 


(A-Zo) 


Substituting  equations  fA-15),  (A-16),  and  (A-20)  into  equation  (A-17)  enables  the 
measurement  error  covanance  matrix  R  to  be  written  as 


R  = 


A,\f  -<,^ 


ft,B,^f  *A,S,q' 


U,B,^f./),8,<r'         Bi^'*8l<!^. 


Lastly,  substituting  the  values  for  A2  and  A3  into  the  above  equation  yields  the  fmal 
expression  for  the  measurement  error  covariance  matrix 


R  = 


g)/    iind^      •*•     CI'    CO 5*" 6, 


co^(e,*e^) 


-vf.  sir^e^co^e,-  <^ ^se,  s,„e^  <r^,  co,\  *  ^^  „^*e, 


cosYe.^e^) 
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APPENDIX  B 
SIMULATION  PROGRAM  LISTING 

V  KALMAN 
[  1  ]  p  THIS  PROGRAM  IS  A  TWO -SENSOR  KALMAN  FILTER  TRACKING  ALGORITHM . 
[2]  fl 
[3]  P 
[L.]  ftTHE  FOLLOWING  IS  A  LIST  OF  THE  PRINCIPAL  VARIABLES  USED : 

[  5  ]  ^APXHAT A  PRIORI  STATE  ESTIMATE  XHAT 

[6]  PiBEl  ,BE2 NORMALLY  DISTRIBUTED  BEARING  ERRORS 

:?]  ^DT TIME  STEP  UN  HOURS) 

18]  ft  ERR ACTUAL  POSITION  ERROR  (TRUE-PREDICTED  ) 

[9]  PC KALMAN  GAIN 

[10]  ftH MEASUREMENT  MATRIX 

[11]  qNTHETAI  ,NTHETA2  .  .NOISY  BEARINGS  FROM  SENSORS  1  AND  2  TO  TARGET 

[12]  P/7Z MEASURED  X,Y  POSITION  (NOISY) 

[13]  PP ERROR  COVARIANCE  MATRIX 

Zl^l  ftPHI STATE  TRANSITION  MATRIX 

[15]  qQ STATE  EXCITATION  COVARIANCE  MATRIX 

[16]  p^ MEASUREMENT  NOISE  COVARIANCE  MATRIX 

[17]  aRNGPCERR PERCENT  RANGE  ERROR 

[18]  ftR2D RADIANS  TO  DEGREES 

[19]  ftTKM TIME  VECTOR  (IN  MINUTES) 

[20]  R  7X  ,  77 TRUE  TARGET  VELOCITIES  IN  X  AND  Y  DIRECTIONS 

[21]  p 7X1 ,772 SENSOR  1  AND  SENSOR  2  VELOCITIES 

.22]  ftXHAT STATE  ESTIMATE 

[  2  3  ]  P  XTK  ,  YTK VECTORS  OF  TRUE  TARGET  POSITIONS 

[2U]  ftXlK»Y2K VECTORS  OF  TRUE  SENSOR  POSITIONS 

[251  R 

[26]  p**************yf*Ayf*****************T'ri*********7"f******yf***i»r** 

[27]  P 

[28]  p  INITIAL  CONDITION  INPUTS: 

[29]  R 

[30]  p  •  INPUT  DESIRED  RANDOM  LINK  I  • 

[31]  i?LSi^7£'^[]/?L^26  506  7 50  0 

[32]  P  ,  ^ 

[33]  '  HOW  MANY  ITERATIONS  ARE  TO  BE  RUN  (K)?* 

[3U]  K^U 

[35]  R 

[36]  •  ENTER  TRUE  TARGET  PARAMETERS :  XTO  .  VXF ,  VXS ,  7^0  ,  77F ,  77S  • 

[3  7]  TRUEST ATE^U 

[3  81  XT Q ^TRUEST ATE  111 

[3  9]  VX^TRUESTATEL2   3] 

[UO]  YTQ^TRUESTATEl^l 

[41]  VY^TRUESTATE  [56] 

[42]  P 

L  4  3  1  '  WHEN  IS  THE  TARGET  GOING  TO  TURN  (WHICH  ITERATION  NO,):^ 


[  u  5  ]  '  NEXT  ENTER  SENSOR  1  AND  SENSOR  2  POSITION  AND  VELOCITY ' 


I] 


[U4]  TURNED 

[U5]  R 

:  -ei 

[U7]  'Z).ir/1  AS  XI  0,7X1, 72  0,772' 

[U8j  SENSPOS-^n 

[U9]  XIO^SENSPOS 

[50]  7X1^S£'/;SP05 

[51]    72  0^5FA75P(9S[3 

[52]    VY2<-SENSP0SL^ 

[53]  R 

[54]  p 

[55]  '  TIME  STEP  TO  BE  USED  (DT.IN  SECONDS )  ?  • 

[55]  DT^n 

[57]  DTS^DT 

[58j  Z)r^Z?rT3500 

[59]  P 

[6  0]  ft 'NO,  OF  POINTS  TO  INCLUDE  IN  THE  REGRESSION:  ' 
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[61]  aNRP^n 

[62]  P 

[63]  •  NOW  ENTER  ALL  A  PRIORI  ESTIMATES  AND  MATRICES :  • 

L5U]  '  ' 

I  6  3  ]  P 

[55]  '  INITIAL  GUESS  XHAT  ( FOUR  ELEMENTS  MUST  BE  ENTERED  )  :  > 

[67"  XHAT-^U 

[5  8]  XHAl^  U  1  oXHAT 

[5  9]  APXHAI^XKAT 

[70]  '  ' 

[71]  p 

[72]  •  INITIAL  P  MATRIX  DIAGONAL  ELEMENTS  (FOUR  ENTRIES  )  :  • 

[7  3]  DIAGONAL^U 

[74]  P^{i\^)o  ,=-^i^)x  i^  n  pDIAGONAL 

[75]  '  ' 

[75]  P 

[77]  a  '  ENTER  OFF-DIAGONAL  POSITIONAL  AND  VELOCITY  COVARIANCES   ( 2  INPUTS  ) : 

[73]  OFFDIAG<-  0  0 

[79]  P[l:3]-?[3;l]-^OFFDi;iC:[l] 

[80]  P[2:4]^PlU;2]-^^FFZ;I/1G[2] 

[81]  PSAVE-^P 

[82]  P 

[33]  '  INITIAL  R  MATRIX  ( u  NUMBERS  IN  THE  ORDER  UL  ,UR ,LL  ,LR) :  ' 

[8U]  /?-f-G 

[8  5]  PS^l^E-e;?^  2  2  pi? 

[35]  p 

[  8  7  ]  p  SET  UP  SOME  VECTORS  NEEDED  TO  GENERATE  THE  Q  MATRIX . 

[38]  RC^((ol).T(Ux360g)I*2 

[89]  rFI-^  12  1  224  2  U 

[SO]  pF1-^F1,-F1 

[31]  rF2'«-  4  4  2  2  4  4  .10p2 

[92]  pF3^(8D4),8tF2 

[93]  RF4f-  43433232 

[94]  rF4'*-F4,F4 

[95]  R 

[95]  '  ' 

[97]  '  ENTER  MAX  BEARING  ERROR  IN  DEC.  FOR  THETAl  AND  THETA2  (^2  INPUTS):' 

[98]  55^0 

[99]  fl 

[101]  R 

R  INITIALIZED  STORAGE  VECTORS  FOR  GRAPHING  PURPOSES  ONLY 

qCVSTORE^OGSTORE^QVERRSTORE^QERRSTORE-^QRESIDERRSTORE'^  ( 2  ,  i^+l )  p  0 
QRNGPCERRSTORE-'rQMSA  VE*-  ( 1 ,  /^+ 1  ;  p  0 

p  THE  FOLLOWING  OUTER  LOOP  IS  FOR  GRAPHING  THREE  SETS  OF  SENSOR  BEARIN 
R  ERRORS  ON  THE  SAME  GRAPH  .  THE  TRACKING  ALGORITHM  IS,  IN  EFFECT ,  BE  IN 
p  RUN  THREE  TIMES  ,  WITH  EACH  RUN  USING  A  DIFFERENT  PAIR  OF  SENSOR 
R  BEARING  ERRORS, 

MAINLOOP^l 

TOP :  BRGERR-^Be  [  ( "  1+2^MAINL00P ) ,  2 ^MAINLOOPl 

P^PSA VE 
R-frRSAVE 

xhat^apxhat 
ai^brgerrli:\ 

A2<rBRGERRL2l 

pSET  up  Q  matrix  initially  to  be  A  4X4  MATRIX  OF  ZEROS . 

g-e-  4  4  po 

R 

Oyt*****A***yr  ******************************************  ********* 

p  THE  PURPOSE  OF  THIS  PROGRAM  SEGMENT  IS  TO 

p  ESTABLISH  TRUE  TARGET  TRACK  (WHICH  INCLUDES  TWO  LEGS )  AND  SET  UP 

p  SENSOR  1  AND  SENSOR  2  TRACKS.  NOISE-FREE  BEARINGS  FROM  EACH  SENSOR 

R  TO  THE  TARGET  ARE  COMPUTED  .  NOISY  ZERO -MEAN  NORMALLY  DISTRIBUTED 

p  BEARING  ERRORS  ARE  GENERATED  AND  THEN  ADDED  TO  THE  NOISE-FREE 

p  BEARINGS  .  FROM  THE  INTERSECTION  OF  THESE  NOISY  LOB  '  S  ,  THE  NOISY 


Ll02 
[103 
[104 
[105 
[105 
[107 
[108 
[109 
[110 
[111 
[112 
[113 
[114 
[115 
[116 
[117 
[118 
[119 
[120 
[121 
[122 
[123 
[124 
[125 
[126 
[127 
[128 
[129 
[130 
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[131]  p  (X ,  y )  POSITION  OF  THE  TARGET  IS  COMPUTED  AND  STORED  FOR  LATER  USE 

[13  2]  ft  IN  THE  TRACKING  FILTER  AS  THE  TARGET  •  S   •  MEASURED  '  POSITION  FOR 

[13  3]  P  THE  KTH  TIME  ITERATION . 

[i3U]  P 

[13  5]  p  SET  UP  SOME  INITIAL  CALCULATIONS : 

[135]  TK^DT^O  ,xK 

[13  7]  TKM-^TKxeo 

[138]  P 

[13  9]  XTKF^XTO+VXZll^DT^O.iTURN 

[140]  XTKS^i    l^XTKF)+VXl22xDTx\iK-TURN) 

[lUl]  XTK^XTKF ,XTKS 

[14  2]  yJAT^Jro  +77  [  1  ]  y^DTx  0  ,  i TURN 

[1U3]  YTKS^i    lfYTKF)+VYL2lxDTx\iK-TURN) 

[144]  YIK^YTKF ,YTKS 

[145]  P 

[146]  xi/^-Hxio+yxixr/: 
[147]  j2A:-ey2  0+7J2xr/: 

[148]  P 

[149]  /?2Z)'e-18  0TOl 

[150]  LVX<rXlK-XTK 

[151]  LVY^Y2K-YTK 

[152]  n 

[153]  THETAlK^{{XlK<XTK)y.ol)^    3o(YTK*  (XIK-XTK+LVX)) 

[154]  r.^Er^2i^-^((y2A:<yr/^)xoi)+  3o(xtk*Iy2K-ytk+lvy)) 

'"155]  P 

[156]  THETAlK<-(LVXxQ.5xol)  +  (^LVX)xTHETAlK 

[157]  THETA2K'*rlLVYx0.5xiol)+(-^LVY)xTHETA2K 

[158]  P 

[15  9]  p  CONTINUE  WITH  CALCULATIONS : 

[160]  P 

[151]  qDATAP0INTS^(2  ,NRP)_pO 

[15  2]  RNGPCERRSTORE^MSAVE^Q 

[15  3]  RESIDERRSTORE-i-NZ*-  2  1  dO 

[164]  GVST0RE^VERRST0RE'^ERRST0REi-BST0RE'^GST0RE^NZFIX<-(2»K-i-l  )pO 

[165]  P 

[166]  X^Ar5r0i?£'-f-4  1  pO 

[15  7]  PSTORE-^K  4  pO 

[168]  p 

[169]  ^-^24pl0000010 

[170]  PHI-e(  i4)o.  =  tH 

[171]  PHI[.1\2^^DT 

[172]  P^J[3;4]-eZ?r 

[173]  ftACCEL*-0.bxDT*2 

[174]  ftGAMMA-^  4  2  pACCEL ,0  »DT ,0  ,0  ,ACCEL »0  »DT 

[175]  J^(i4)o.rii; 

[176]  P 

[177]  p ************* ******************************************* 

[178]  P 

[179]  p  START  MAIN  LOOP: 

[180]  P 

[181]  N^l 

[182]  P 

[183]  P 

[184]  LOOP :  (;-^P+  .  X  ( <5ff )  +  .  X  S/?+f/+  .  xp+  .  X  J5ff 

[185]  P 

[185]  P 

[187]  p  'KALMAN  GAIN:  ' 

[188]  P 

[189]  GVST0RELl',N']frGl2:i: 

[190]  GVSTOREL2',Nl'f-Gl^;2: 

[191]  GSTOREll',N'j*-GLl;l'] 

[192]  GSr(?/?£'[2;/\/]-eG[3;2] 

[193]  PC 

[194]  R 

[195]  R 

[196]  P-«-(J-C+.x^)+.xP 

[19:^]  PSTORE^PSTORE ,L0  .b^l  +  KNl  P 

[193]  R 

[199]  P 

[200]  p  '  UPDATED  P:  ' 
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[201]  rP 

[202]    ZK^  2  1  ^XTKIN-}  ,YTKIN'] 
[203]  UNTRUE  TARGET  POSITION^ 
^204]  pZA' 

[ 2 05  J  fi  THE  NEXT  SECTION  GENERATES  2  NORMAL  BEARING  ERRORS . 

[207]  R 

[208]    £/-f-(?2clO*10  )^10*1C 

[209]    SSQ^i    2^®:/[l]  )*r     - 


[210]    S?R1^A1t(R2D^1.^6) 
[211  j    SPR2^A2rlP.2D^1.96) 


[212]    E  STORE  [  1 ;  A']  ^HEl  ^5P/?1  ^SSQ  x2o(^[2]x2xol) 

[213]    5SrO/?£ [2  ;iV]^5£:2^5Pi?2xSS$xio  (a [2 ]x2xol) 

[214]   fl 

[215]    NTHETAl-^rTHETAlKZNl-hBEl 

[215]    NTHETA2'<rTHETA2KLNj+BE2 

[217]  R 

[216]    S^  1-^1  o /vry ETA  1 

[219]    CAl*-2oNTHETAl 

[2  20]    SA2^1oNTHETA2 

[221]    CA2'^2oNTHETA2 

[22  2]    D<r2oiNTHETAl+NTHETA2) 

[223]   R 

[2  24]    /VZFJX[l;W]-^;VZ[l;]^C(y2A:[W]xSA2xCAl)-(XlA[/7]xS>^lxS/12})TD 

[22  5]    A'2FIX[2;A7]'f-/VZ[2:]-^((XlA[A/]xSAlxCi42)-(y2A[A']x5A2xSAl))TD 

[226]  f^NOISE LliNl ^XTK INl  -NZFIX  [  1 ;  A'] 

[227]  aNOISEL2;N:i^YTK\:N2-NZFIXL2iNl 

[228]  fl  '  TRUE  TGT  X  yY  POSIT  MINUS  NOISY  X  ,Y  POSIT' 

[229]  R  NOISE L;N1 

[230]  R 

[231]  R 

[232]  fl 

[233]    RESIDERR-^NZ'H+.yXHAT 

[234]    XHAT^XHAT+G+.xRESIDERR 

[235]   R 

[235]  R 

[23  7]    XHATSTORE-^XHATSTORE ,XHAT 

[238]  R 

[239]  RESIDERRSTORE-^RESIDERRSTORE ,  RESIDERR 

[240]    R 

[241]  uDATAPOINTS^DATAPOINTS ,NZ 

[242]  f^DATAPOINTS-^DATAPOINTS ,  10  10  /Xf/AT 

[243]   R 

[244]  R DATAUSE<r  ( 2  ,  -A/PP )^DATAPOINTS 

[245]  RSJWX^  +  /MrAa5£'[l;] 

[246]  RSy«y'^  +  /Z?ArA(/5E[2;] 

[247]  RSaWX2'^  +  /DArA/7S£:Ll:]*2 

[248]  RSa?/:xy-^+/I5Ar>i£/SE[l ;  1  xCAr>!iy5E[2  ;  ] 

[249]  QDENOMINATOR'<riNRPySUMX2)-SUMX*2 

[2  50]  qM^((NRP^SUMXY)-SUMX^SUMY)*DENOMINATOR 

L2512  ^MDEG^R2Dx(ol)+    3oM 

[2  52]  RM5A7E-eM5/17E,WZ)EG: 

[253]  fl 

[254]    7E/?PSrOP£:[l;A7]-«-yX[l+A7>r^PA/+l]-XM2'[2;] 

[255]    VERRSTORE L2;N1'^VYl1  +N>TURN+1j  -XHATl^;  ] 

[255]    RTRUE'^i+7^ZKxZK)^0.5 

[2  5  7]    PX-f-  10  10  /XHAT 

[253]    P//AI'<e(+/PXxPX)*0.5 

[259]    RNGPCERR^IQ 0 x ( tRTRUE ) x \RTRUE-REAT 

[250]    RNGPCERRSTORE^RNGPCERRSTORE , RNGPCERR 

[251]  R 

[252]  fl 

[253]  R  '  UPDATED  XHAT:  ' 

[254]  fl  XHAT 

[255]  EPP^Z/^-  2  1  D  ,XffAr[l  3  ;] 

[255]  fl  '  ERROR  VECTOR:  ' 

[257]    PPPS2'CPP[l;A'j-^PPP[l;] 

[258]    PPPSrOPE[2;A7j^EPP[2;] 

[259]  f^ERR 

[270]  fl 
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C271] 

[272] 

[273] 

[274] 

[275] 

[276] 

[277] 

[278] 

[279] 

[280] 

[281] 

[282: 

[283: 

[284" 

[285] 

[286" 

[287 

[288 

[289] 

[290: 

[291' 

[292 

[293] 

[294] 

[295] 

[295" 

[297 

[298 

[299 

[300 

[301 

[302 

[303] 

[304: 

[305" 

[306] 

[307" 

[308 

[309 

[310 

[311 

[312 

[313] 

[314: 

[315] 

[316] 

[317] 

[318] 

[319] 

[320] 

[321] 

[322] 

[323] 

[32U] 

[325] 

[326] 

[327] 


■Q-^PHI+,xp-\-,xtSiPHI 


' PREDICTED  P : 
P 


Q 
P 
P 
fl 

n 

A 

XHAT<-PHI+.^XHAT 
fl 

R  '  PREDICTED  XUAT:  ' 
fi  XHAT 

QXl-frXHATLl',! 
RXn<-XHATl^;2 
flX^.-s-X4,X4,X2,X2,X4,X4,X2,X2,X2,X2,X2,X2,X2,X2,X2,X2 

pXB^X4,X4,X4,X4,X4,X4,X4,X4,X4,X4,X2,X2,X4,X4,X2,X2 

(^QTERMS'<rC^FlxXA>iXB^DTS*F^ 

oC-f-  4  4  oQTERMS 

R  CALCULATE  TERMS  IN  THE  R  MATRIX : 

/?l<-(X;:/;i2'[3;]*2)  +  (X^Ar[l:]-Xlii:[A^]  )*2 
5IG12-«-/?lx(3o.^lfi?2Z?xi.96  )*2 
R2'<riXHATLl;:i*2)+iXHATL3:J-Y2KLNl  )*2 
SJC2  2-c-/?2x(3oA2*i?2Dxl.96)*2 

R 

Dl'(-D*2 

ULf-((SIC12^SA2*2  )+(SJC22xCAl*2 

L/?^((SJC12xCA2*2)  +  (5I(;22xS/41*2.  . 

CROSS-^-i  iSIG12^SA2^CA2)  +  (SIC22xSAlxCAl))*Dl 

R^2  2  pUL  .CROSS , CROSS »LR 

R  •  R  MATRIX : ' 

fl  R 

R  '  » 

R  •  ITERATION  NO,   '  ,^iN+l) 

R  '  ' 

N^N+1 

^{N<K+1)/L00P 

R 

R ****************************************************** 

LQM<rO,5^1  +  l<MAINLOOP 

' GVSTORE-'rQGV STORE  .ILQMI  GVSTORE 

GSTORE^QGSTORE .lLQMI  GSTORE 

VERRSTORE^rOVERRSTORE .  \.LQM'\  VERRSTORE 

ERRSTORE^rQERRSTORE  .LLQMl  ERRSTORE 

RESIDERRSTORE-^QRESIDERRSTORE.lLQMI  RESIDERRSTORE 
..NGPCERRSTORE^ ( 1 .oRNGPCERRSTORE )oRNGPCERRSTORE 
QRNGPCERRSTORE-f-QRNGPCERRSTORE,  LLQMl  RNGPCERRSTORE 
fi  MSA  VE"-  ( 1  .  p  MSA  VE)  p  MSA  VE 
uQMSAVE^QMSAVE ,  ILQM']  MSAVE 
MAINLOOP^MAINLOOP-^1 
URL-rRLSAVE 

■*iMAINLOOP<>Q  ,  5xpB6  )/TOP 
'    COMPLETE, ' 


))*D1 
))*Z?1 
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